{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import sklearn\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.model_selection import train_test_split"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "class LinearRegression:\n",
    "    def __init__(self, learning_rate=0.0001, n_iters=30000):\n",
    "        self.lr = learning_rate\n",
    "        self.n_iters = n_iters\n",
    "        self.weights = None\n",
    "        self.bias = None\n",
    "\n",
    "    def fit(self, X, y):\n",
    "        n_samples, n_features = X.shape\n",
    "\n",
    "        # init parameters\n",
    "        self.weights = np.zeros(n_features)\n",
    "        self.bias = 0\n",
    "\n",
    "        # gradient descent\n",
    "        for _ in range(self.n_iters):\n",
    "            # approximate y with linear combination of weights and x, plus bias\n",
    "            y_predicted = np.dot(X, self.weights) + self.bias\n",
    "\n",
    "            # compute gradients\n",
    "            dw = (1 / n_samples) * np.dot(X.T, (y_predicted - y))\n",
    "            db = (1 / n_samples) * np.sum(y_predicted - y)\n",
    "            # update parameters\n",
    "            self.weights -= self.lr * dw\n",
    "            self.bias -= self.lr * db\n",
    "\n",
    "    def predict(self, X):\n",
    "        y_predicted = np.dot(X, self.weights) + self.bias\n",
    "        return y_predicted\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "class LinearRegressionWithSGD:\n",
    "    def __init__(self, learning_rate=0.0001, n_iters=30000):\n",
    "        self.lr = learning_rate\n",
    "        self.n_iters = n_iters\n",
    "        self.weights = None\n",
    "        self.bias = None\n",
    "\n",
    "    def fit(self, X, y):\n",
    "        n_samples, n_features = X.shape\n",
    "\n",
    "        # init parameters\n",
    "        self.weights = np.zeros(n_features)\n",
    "        self.bias = 0\n",
    "\n",
    "        # stochastic gradient descent\n",
    "        for _ in range(self.n_iters, batch_size=5):\n",
    "            # approximate y with linear combination of weights and x, plus bias\n",
    "            y_predicted = np.dot(X, self.weights) + self.bias\n",
    "            \n",
    "            indexes = np.random.randint(0, len(X), batch_size) # random sample\n",
    "        \n",
    "            Xs = np.take(X, indexes)\n",
    "            ys = np.take(y, indexes)\n",
    "            y_predicted_s = np.take(y_predicted, indexes)\n",
    "            \n",
    "            # compute gradients\n",
    "            dw = (1 / batch_size) * np.dot(Xs.T, (y_predicted_s - ys))\n",
    "            db = (1 / batch_size) * np.sum(y_predicted_s - ys)\n",
    "            # update parameters\n",
    "            self.weights -= self.lr * dw\n",
    "            self.bias -= self.lr * db\n",
    "\n",
    "    def predict(self, X):\n",
    "        y_predicted = np.dot(X, self.weights) + self.bias\n",
    "        return y_predicted\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'lr': 0.0001, 'n_iters': 30000, 'weights': array([ 0.53424952,  0.33346318, -0.00779044,  0.09959714,  0.2180196 ,\n",
      "        0.06196114,  0.13746151,  0.00133901]), 'bias': 0.03227018649613295}\n",
      "0    -1.288078\n",
      "1    -0.919539\n",
      "2    -1.007862\n",
      "3    -0.795213\n",
      "4    -1.324515\n",
      "        ...   \n",
      "92    0.873505\n",
      "93    0.423971\n",
      "94    1.621920\n",
      "95    1.637418\n",
      "96    1.469376\n",
      "Name: lpsa, Length: 97, dtype: float64\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWoAAAD4CAYAAADFAawfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAmIElEQVR4nO3dd3hUZdo/8O+TZEJCDSW0QAhFQ4dABAQsuKuooKKIFVd/i4uIYl0U1t1V13cRRXfXjqztXUIXRFAW9BUsIIIJCT0gHUJJKKGlz9y/P0hiSM7MnJnMzCnz/VyXl5pMZu6TwHee3OcpSkRARETmFWF0AURE5BmDmojI5BjUREQmx6AmIjI5BjURkclFBeNJmzVrJklJScF4aiIiW8rIyDguIvFanwtKUCclJSE9PT0YT01EZEtKqf3uPsfWBxGRyTGoiYhMjkFNRGRyDGoiIpNjUBMRmZyuWR9KqX0AzgJwAigTkdRgFkVERL/yZXreEBE5HrRKiIhIE1sfRGQrxWVOTFm2HYfzC40uJWD0BrUA+EoplaGUGqv1AKXUWKVUulIqPS8vL3AVEhHpdDi/EHe8/xNmfL8Hq3bkGl1OwOhtfQwWkRylVHMAXyulskXk+6oPEJEZAGYAQGpqKk8jIKKQWrPrOCbMyURJmQvTR/fF9d1bGl1SwOgaUYtITvm/cwF8BqBfMIsiItJLRPDut7tw34fr0LReND5/dJCtQhrQMaJWStUDECEiZ8v/+zoAfwt6ZUREXpwpKsUf52/EV9uOYXjPVnhlZE/UqxOULYwMpeeKWgD4TClV8fjZIrI8qFUREXmx4+hZjEvLwMGTBfjr8K74f4OSUJ5TtuM1qEVkD4BeIaiFiEiXz7NyMGnhZtSPicLsPwxAv/ZNjC4pqOz3OwIR2Vap04Upy7bj4zX7cFlSY7xzTx80bxhjdFlBx6AmIkvIPVOE8bM2IH3/Kfx+UHtMvrEzHJHhsRSEQU1Eprd+70k8MnsDzhWV4c27U3Bzr9ZGlxRSDGoiMi0RwUdr9mHKsu1IbFIXaWP6I7llA6PLCjkGNRGZ0vniMjy7cBO+2HQE13Vtgdfu6IWGMQ6jyzIEg5qITGd33jmMm5mB3Xnn8Oz1nTHuqg62nXqnB4OaiExl+Zaj+OOCjYiOisDMMf0xqFMzo0syHIOaiEyhzOnCa1/txPTvdqNX2zi8d28ftI6LNbosU2BQE5Hhjp8rxmNzMvHj7hO4t38i/npTV9SJijS6LNNgUBORoTIPnML4WRtw8nwJpt3eE6NS2xpdkukwqInIECKCWesO4MWlW9GyUQwWPjwQ3RMaGV2WKTGoiSjkCkuceG7xZizakIOrk+Pxrzt7I65utNFlmRaDmohC6sCJAjyUloHso2fwxG8vwWPXXIKIiPCdeqcHg5qIQmZVdi4en5sJAPjo/sswpHNzgyuyBgY1EQWd0yV445tf8OY3v6Brq4aYProvEpvWNbosy2BQE1FQ5ReU4PG5WfhuZx5G9mmDv9/aHTEOTr3zBYOaiIJmS85pjEvLQO6ZYvz91u64p19iWC8F9xeDmoiCYn76Qfxl8RY0qReN+eMuR++2cUaXZFkMaiIKqOIyJ15Ysg1z1h/AwI5N8dbdKWhav47RZVkag5qIAiYnvxDj0zKw8dBpPHx1Rzx97aWICpNTWIKJQU1EAbH6l+N4bG4mSstceP++vhjaraVfz7M4MwfTVuzA4fxCtI6LxcShyRiRkhDgaq2FQU1EteJyCd77bjde/2oHOjWvj+mj+6JDfH2/nmtxZg4mL9qMwlIngAsj9MmLNgNAWIc1fychIr+dKSrFQ2kZmLZiB4b1bI3Pxg/yO6QBYNqKHZUhXaGw1IlpK3bUtlRL44iaiPyy4+hZPDQzHYdOFeL5m7rigYFJtZ56dzi/0KePhwsGNRH57POsHExauBn1Y6IwZ+wAXJbUJCDP2zouFjkaoRzuBwgwqCks8AZVYJSUuTBl2XZ88uM+XJbUGO/c0wfNG8YE7PknDk2+qEcNALGOSEwcmhyw17AiBjXZHm9QBcaxM0V4ZNYGpO8/hTGD22PSDZ3hCPDUu4qfh9XeVIM9EGBQk+15ukFl9gAwi3V7TuCR2ZkoKCnDW3en4KZerYP2WiNSEiz1cwnFQICzPsj2eIPKfyKCD37Yg3s+WIeGMVFY/MigoIa0FYVipgpH1GR7vEHln/PFZXhm4SZ8uekIhnZrgddG9UKDGIfRZZlOKAYCHFGT7U0cmozYattq8gaVZ7vzzuGWd9bgv5uPYNINnTF9dF+GtBvu3vADORDQHdRKqUilVKZS6ouAvTpRCIxIScDLt/VAQlwsFICEuFi8fFsPS/VBQ2n5liO45e01OHW+BGlj+mPcVR25NakHoRgI+NL6eBzAdgANA/bqRCFitRtURihzujDtqx14/7s96NU2Du/d24ftIR1CMVNFV1ArpdoAGAbg7wCeCtirE5EpHD9XjAmzM7F2zwmMHpCIvwzvijpRPIVFr2APBPSOqP8F4BkADdw9QCk1FsBYAEhMTKx1YUQUGhsOnML4tA04VVCC10b1wu192xhdElXjtUetlBoOIFdEMjw9TkRmiEiqiKTGx8cHrEAiCg4Rwcyf9uPO99fCEaWwaPxAhrRJ6RlRDwJws1LqRgAxABoqpdJEZHRwSyOiYCksceK5xZuxaEMOurZqiFPnSzD8zdWWWQkYbrwGtYhMBjAZAJRSVwP4I0OayLr2nziPh2ZmYMexs7i+e0t8m52LojIXAC6v18OIfWM4j5oojHyz/RiGv7UaR04X4aMHLsPmQ6crQ7oC9392r2K5eE5+IQS/vrEtzswJ6uv6tDJRRL4F8G1QKiEir/wdzTldgje++QVvfvMLurVuiOmj+6Jtk7pcXu8jo/aN4RJyMjVuT/orfzf/OXW+BI/Py8L3O/Nwe982+J8R3RFTvkCDy+t9Y9QbG1sfZFpG/ZppVu5Gc0/My8KgqSs1vy+bD53G8LdW46fdJzDl1h6YdnvPypAGuLzeV6FYLq6FQU2mxfPzLuZp1Kb1Jjb/54MYOf1HiAgWjLsc9/RPrLEUnMvrfWPUGxtbH2Ra7J9ezF2bokLFm9j13VvixaVbMWf9QQzq1BRv3pWCpvXruP06Lq/Xz6iDDRjUZFrh0j/V24fXOqaqupz8Qtzx/lpsOnQa46/uiKevS0ZkBDdUCiQj3tjY+iDTCof+qS99+KptCnciFLA37zxm3NcXz1zfmSFtEwxqMq1w6J/62ocfkZKANZOuwb/u7F3jTQwAmjeIweePDsJ13VoGpV4yBlsfZGp275/624ev+J688t9sHDlTBADokxiHtAf7o240/1rbDX+iRCFWtScdoRScIjUeo6cP37lVA9RxRCAqQuG5YV3wwMAkbvBvUwxqohCqvmhFK6T19OE/z8rBsws3oWGMA3PHDkBqUpOg1EvmwKAm8iKQqyO1etIAEKkUXCJen7+kzIUpy7bjkx/3oV9SE7x9bwqaN4jxq5ZgMuuKUrPW5Q2DmsgDf5dtu+Ou9+wSwd6pwzx+7bEzRRg/awMy9p/CmMHtMemGznBEmm8+QKC/Z3avSw/z/ZSJTOTFpVsDujrS3yXIP+05gWFvrsb2I2fw9j0p+MvwrqYMacC8K0rNWpce5vxJE5nA4swcnCoo1fycv6sjfZ0bLiL44Ic9uPeDdWgYG4XPHxmE4T1b+/XaoWLWFaVmrUsPtj6I3PA00vJ3daQvS5Dnrj+AF5duQ2GpEzGOCIy9ogMuaeH22FLTMOuKUrPWpQeDmsgNTyOt2qyO1DM3/P3vdmPqf7NRMSekqNSFF5duQ4wjMuT9VF9vwGktdTfDilKz1qUHWx9EbrgbacXFOoIali8u2YqXq4R0BSP6qf5sNWvWFaVmrUsPJRrzOGsrNTVV0tPTA/68RKH058WbkfbTgRofHz0gEf8zokfAX6/M6cLYmRlYmZ3r8XEJcbEhm142aOpKzXZBQlws1ky6JmivG46UUhkikqr1ObY+iNxYlZ3n08drM0c372wx7v73T9iVe87rYyuCMxTTy6x8A85OGNQU1jyFqy8hVZs5uhn7T+H3n/yM04XaM0w8CfZ5fVa+AWcn7FFT2PLWf/VlzrM/c3RFBDPX7sNdM9biXHGZ39cRzNFtOGw1awUMagpb3sLVl5DytUVQWOLE0/M34i+fb8UVl8TD6fL/XlEwR7dWvgFnJ2x9UNjyFq6+zHn2pUWw/8R5PDQzAzuOncVT116KR4d0whWvrtL8+kilEOOIwPkS7VNdAj26ddcKYjAbi0FNYUtPuFYNqYoQe3JeFhrFOqAUkF9QitZxsRjSOR4LM3K8ztH9v23H8OT8LEQohY8fuAxXJzcH4H6O78i+CZi3/qBm/Y3rOvD8Td0CFqJW3gvD7tj6oLA1pHO87o9X72fnF5biVEFpZW97YUYORvZNcNsicLoEr3+1Aw/+Jx3tmtbFFxMGV4Y04L7FsCo7D6UabZG4WAcy/3pdQAPUynth2B1H1BS2fJl+52570gqFpU6sys7TnFt86nwJHp+Xhe935mFU3zZ4aUR3xGgco6XVYnhyXpbm6/kzQ8QbTsUzLwY1hS1fgklPWGk9ZvOh0xiXloG8s8V4+bYeuOuytj6dwhLK6XGcimdebH1Q2PJl+p2esKr+mHk/H8DI6T9CRLBg3OW4u1+iz0dlhXJ6HKfimReDmnyyODMHg6auRPtJX2LQ1JUe93wwO1+CSeux7r6uqNSJyYs24dmFm9EvqQm+eOwK9Gob51eNoZwex6l45sW9Pki36rMCgAsBZeW/zL4s+6762OqzPiq+7tCpAjyctgGbc07jkSEd8dS1yYiM4IGz5J2nvT68BrVSKgbA9wDq4EJP+1MRed7T1zCo7Ykb9Hj2/c48PDY3E06n4PU7euG6bi2NLokspLabMhUDuEZEzimlHABWK6X+KyI/BbRKMj3OCtDmcgne/XYXXv96Jy5t3gDT7+uL9s3qGV0W2YjXoJYLQ+6KLb0c5f8Evl9CpsdZATXNXrcff1u6DUVlLsQ6IvH7wUkMaQo4XTcTlVKRSqksALkAvhaRdRqPGauUSldKpeflac9PJWvjrICLvbtqF577bAuKylwALsylfmHJNkveYLXTTWI70hXUIuIUkd4A2gDop5TqrvGYGSKSKiKp8fHaK77I2jgr4FefZR7CtBU7THEKS235c4oLhZZPC15EJF8ptQrA9QC2BKckMrNw36CnpMyFv3+5Df+7dr/bx1itZ+9p6Xg4/6zNxGtQK6XiAZSWh3QsgGsBvBL0yohM5ujpIoyflYENB/Jx9aXx+H5nHlwajwt2z742J8lo4U1i89Mzom4F4H+VUpG40CqZLyJfBLcsInNZu/sEJszZgIISJ+6/vB3mpx/SDOlg9+yDscMdbxKbn55ZH5sApISgFrKQQI/qzEpE8MEPezF1eTbaNa2LOX8YgAc+/llzg6ZIpYLes3fXpnhhyVa/fx7utlgN15vEZsRNmchndtm32NubzbniMjzz6UYs23wUN3RviVdv74kGMQ63LQGXSNCv391r5xeWIr98Rz1ffx6+HJBAxmBQk8/scPPJ25vNrtyzeGhmBvYeP4/JN3TG2Cs7VG6oZGSrwN1rV+frzyPcbxKbHYOafKb35pOZ2yPu3myenr8RWQfzsSD9IGIckUh7sD8Gdmx20eOMbBVovbY7vBloHwxq8pmeEaUv7REjAt1diDlF8MmP+9CuaV3MHTsArRrFatY4sm8CVmXnhfxNSKtNUVBShlMFNQ8SaBTrCHo9FBoMaqqkNzD1jCj1tkdq2+/2N+S9tRBKy1wXhXT1Ghdm5Bi22Kd6m2JxZg4mLthY48iu8yVlWJyZY5rfYsh/3I+aAPi2Ok3PCkW97ZHanNNXmxV13vaXPnK6KCA1hsKIlATUj6k55ip1imlqpNrhiJoA+H6D0NvNJ7033Gqz2MKXmrVG3lNu7Y6n5m/U3GGsap1WWBCSr9H6AMxVI/mPI2oCEPgw0ruBk7uZEhFKed0gyJebmtVH3pMWbsKsdQcgAKrv61+9Tl+O7DKKFWok/zGoCUDg/6Lr3cDJXQvCKeK1neGt5ood4Z6Yl1Vj5F1U5kL6/lN4+tpL8dqoXh7rtMKugVaokfzHo7gIgLHHbFVtS0QoBafGn0mtU2Q81QxA1zS2fVOH+Vyj2aYaVrBCjeRerY7i8geD2prM8Be9/aQvNXvGCsBejVB1V7O7Y8Oq4hFiZCa1PYqLwoQZVqf5uurPXc3eeuuOSBXQtoAZ3uTIvtijJlMJVK/VW2+9XnSUxyD15cQTbrxPwcagJlMJ1CkyE4cmwxGp3H7+dKH2dDZAO3ifnJeFJDehbfZ51mR9bH2Q6dS2BVNU6sTa3SdQ6nR//8XTiFsreCueSWvlpBXmWZO1cURNtnLoVAFGTV+LeekH8eiQTvjHqF6a0/8KypdXa/EWsNVHy6GYw8zDZ8MbR9RkG9/tzMPjczPhdAr+/btUXNu1BQAgIkLhhSVbK/drBoBTBaVu9xTRs5Vo1TAP9m56dtn/m/zHETVZnssleOubX/DAx+vRsmEMlkwYXBnSwIUwq1en5pjEXR/Z2z4gwMWj5WCfzs4eOHFETZZ2urAUT83LwjfZuRjRuzWm3NYDdaNr/rH2pY9cdSvRnPxCKOCiud1ao+VgTm1kD5wY1GRZ24+cwbi0DOScKsSLN3fD7y5vV3kKS3W1mZ9t9BxpHj5LDGqypM8yD2Hyos1oFOvAvIcGoG+7Jh4fX5s+stELgXj4LDGobc7o0WCgfZp+EM8v2YrzJU5ER0bgkas7eQ1p4ELYpu8/iTnrDsIpgkilMLJvYAI42N9jHj5L3OvDxozcaCkYPlmzFy9+sQ1V/8jqvR6t70VVjes68PxN3Xz+vtjte0zG8bTXB2d92JhRswWCMed37e4T+Fu1kAb0X4/W96KqUwWlmPjpxotq1XMdnJFBocDWh415mi0QrF/XAz3nV0Tw7x/24JXlO+By88ufntkPeh5TcXTViJQE3dfBGRkUChxR25i7WQGNYh1B20QokCPMc8VlGD9rA6Ysy8bQbi3QqlGM5uP0zH7QO0OiImD1XgdPVqFQYFDbmLud6JRC0H5dD9QIc1fuWdzy9mqs2HoUf7qxM965pw+evb5zjetRAIZ0jvf6fHoWsQC/Bqy7lYnVr4Mnq1AoMKhtzN2KuUAfhFq1lxvhYR6zXl9uOoJb3l6D04WlSHuwP8Ze2RFKKYxIScDIvgmo+goCYGFGjtffBqp+L9yp2KN6cWYO3O27V/06gr0qkQhgj9r2tOYAV6y4q86fX9er93K1jtHSO8Isc7rwyvJs/PuHveiTGId37+2LllXaHYszcyoPpK3K02npVVVfxFJ1/4+qsz4GTV3p9pQZreswep412R+DOswszszB+eKyGh/399d1d7MpIpWCS0TzRqXWjcxBnZrh0dkbsG7vSdx/eTs8N6wroqMiLvqaiQs21pj1UcHX3wY8hau7toeAmyCRMRjUJhLshRPu5hL7O4cYcB+QLpHKMw4rWiOH8wvRKNaB8yVllXtF5+QX4plPNyE2OhLFZU78885euDWlTY3nm7ZiB0rdTftA4G7eVbQ9tF7JU9uEKJi8BrVSqi2A/wBogQt/fmeIyBvBLizchGIrS3ej37pejqXyxNs+FNWvK1/jZJUSpwvOYsEXEwajS6uGmq/jbcQcqJt301bs8KntQRQKem4mlgF4WkS6AhgA4BGlVNfglhV+QrFwwl3Y5eQX+r1AxdusB28LTSo4XeI2pAHPI+a4WEfA3szY9iAz8hrUInJERDaU//dZANsB8E9sgIVi4YS7sFOA33Oqvc160Fu/t7bCxKHJcETUnIvhiFR44eZuul5Dj0g3s1bcfZwoFHyanqeUSgKQAmCdxufGKqXSlVLpeXl5ASovfIRi4YTW6FerH+tuJO9uSfWIlASsmXQN/nlnbwDAk/OyKj+vp349NzJHpCRg2qheiIt1VH6scV0Hpt3eK6AjXa1ZK54+ThQKum8mKqXqA1gI4AkROVP98yIyA8AM4MKmTAGr0Kaq3zgc0jkeCzNygrqVpdYubHoXdnjrobv7/Mi+CTWuSwFoEBOFs0VlPt00DcU0uAQ33xPeSCQj6QpqpZQDF0J6logsCm5J9qcVagszcjCybwJWZecFdSvL6mE3aOpKXXOqPfXQR6QkuP38quw8PDesC176YhuKy1yoGx2JF2/uhlGpbQN4VYHDvZ/JjPTM+lAAPgSwXUT+EfyS7M9TqK2ZdE1Ia9EbTN566J5uVL737W4IgKm39cBd/RIDU3iQcO9nMiM9I+pBAO4DsFkplVX+sT+JyLKgVWVzZtpxTW8wuWuTRChVuXTcUx/303GXo2ebuIDWHixcaUhm4zWoRWQ14HbrA/KD2c7A0xNMWiNv4NebbO5COrllA8z5wwA0qRdd+TG7nTpDFGzclMkAVttxrSJYC0udldPU3E1Xq/oH6rquLbDssStqhHSwtlglsisuITeAlfqgWpsuxToi3S5iceHCjI5/3tEbv+3aosbnvd2UJKKaGNQGsUof1F2wRrrpSUdFKCx9dDCSmtXTfD4z9eeJrIJBrVO49lXdBajWyDpSKUy5tYfbkAbM158nsgL2qHUI576quwBNiIvFo0M6IbJ8WXejWAdeG9UT0VERHvcNsVp/nsgMGNQ6hPNJ0+6C9erkeLy16hc0qx+NhQ8PxMbnr4NSyusbGk9EIfIdWx86hHNftfqNz1aNYtAhvj5mrTuA/u2b4O17+iC+QZ3Kx+i5UWiV/jyRWTCodQj3vmpFsB45XYjxszZg9a7jGHtlBzwzNBlRkb/+UhbOb2hEwcTWhw7sqwI/7j6O4W+uxs6jZ/HuvX3wpxu7XBTSQGh2ACQKRwxqHcK5ryoieP+73Rj9wTrE1XXg80cH4cYerTQfyzc0ouBg60Mnu/dVtaYf/qZLc0xcsAnLtx7FjT1a4tXbe6F+Hfd/ZKy0kIfISpQEYUP01NRUSU9PD/jzUnBoHXpbJyoCjWIdOHG+BJNv6Iwxg9tD8ZQToqBRSmWISKrW5ziitjk9C3W0ZmsUl7lw/FwxZv9hAAZ0aBrKkomoGga1jek92dzdrAyXgCFNZAK8mWhjehfqeFp9SETGY1DbmN55zSP71LzZx9kaRObBoLYxb/OaRQSfrNmLd7/djWb1o9G8QZ2wm35IZAXsUVuQ3p38PJ2HWFBShsmLNuPzrMP4bZfmeP2O3mgU6wjlZXhl5h0LzVwb2Q+D2mL03iCs+v/VA6VX2zjc+s6P2Jl7FhOHJuPhqzoiIiK0U++8BZ0v1xlqZq6N7InzqC1m0NSVmvuOJMTF6jrB/KutR/H0/I2IilR4464UXHlpfDDK9Ehr3nasI/KidkttrzOYzFwbWZenedTsUVuMvxsfOV2CaSuyMXZmBtrH18PSCYMNCWlA32wUM2/wZObayJ4Y1Bbjz8ZHJ8+X4P6P1uOdVbtxd7+2mP/Q5WjTuG6wSvTKXaBVHaWaeYMnM9dG9sSgthhfNz7aeDAfw9/8Aev3ncQrI3vg5dt6Iqba14eau0BTQOUhA2be4MnMtZE9MagtRu9OfiKCOesPYNT0tVBKYeG4gbjzskRjiq5m4tBkaN26FKCy/WHmHQvNXBvZE28m2lBRqRN//XwL5qcfwpWXxuONO3ujcb1oo8u6SNKkLzU/rgDsnTostMUQmQA3ZQojB08WYFxaBrYePoPHrumEx397aeUBtGaSEOan5hD5gkFtI9/uyMXjc7PgEsGH96fiN11a+P1cwV7Q4WkxDhFdjEFtAy6X4O1Vu/DP/9uJ5BYN8P59fdGuaT2/ny8UCzp4yACRfgxqiztdUIon52dhZXYubk1JwJRbeyA2unazOrydJh6o0bbdT80hChQGtQVVBGVOfmFl//mlW7ph9IB2ATmFxdOCDi6fJgo9r9PzlFIfKaVylVJbQlEQeVYRlBU34pwuQaRSaBDjCNhRWZ4WdOjd45qIAkfPPOpPAFwf5DpIp1eXZ9cIyhKnK6BB6WlBB5dPE4We16AWke8BnAxBLeTFkdOFOHy6SPNzgQxKTws6uHyaKPQC1qNWSo0FMBYAEhPNsQLOTn7cdRwT5mRC4cIKvuoCHZTubvRxWh1R6AVsCbmIzBCRVBFJjY83Zlc2OxIRTP9uN0Z/uA6N60Vj0g2dDd1ngsuniUKPsz5M7GxRKSYu2ITlW49iWI9WeOX2nqhfJwotGsbonh4XjIUrnFZHFFqmCWp/AsXOxyHtPHYW42ZmYP/JAvx5WBeMGdy+claH3qDkVDoie9AzPW8OgLUAkpVSh5RSYwJdRNUpZ4JfA6Viy8tAfY1VLN14GCPeWYMzRWWY/WB/PHhFB7+m3nEqHZE96Jn1cbeItBIRh4i0EZEPA12EP4FixxAqdbrw0hfbMGFOJrq0aogvHxuM/h2a+v18nEpHZA+maH34Eyh2C6Hcs0V4dFYm1u87iQcGJuFPN3ZBdFTt7vW25g51RLZgiqD2J1CsHkJV++tN60ejpMyFUqfgjbt645begekfcyodkT2Y4oQXf442GtJZewqgu4+bSfX++vFzJThbVIYJ13QKWEgDnEpHZBemGFH7s+Xlquw8nz5uJlr9dQEwa90BjB/SKaCvxal0RNZniqAGfA8UK/eotVo2gDVqJ6LQM0Xrwx9W3XPiq61HNQ92BcxfOxEZw7JBPXFoMhzVzgJ0RCjT3ihzugSvLs/G2JkZaNM4FnWqzejgTT4icseyQQ0ANYam5jvDFQBw4lwx7v9oPd79djfu7peIr5+6Cq+M7MmbfESki2l61L6atmIHSp0X7yNX6pTK46LMIutgPsanZeD4+RK8OrIn7risLQDe5CMi/Swb1Ga/mSgimLP+IF5YshXNG9bBoocHontCI6PLIiILsmxQm3nBS1GpE39ZvAULMg7hykvj8cadvdG4XrTRZRGRRVm2R+3PIplQOHiyACPf+xELMg7hsd9cgo8fuIwhTUS1YtkRtT+LZIJt1Y5cPDE3CyKCD+9PxW+6tDCsFiKyD8sGNWCeG3Iul+DNlb/gjW9+QeeWDTF9dB+0a1rP6LKIyCYsHdRmcLqgFE/My8SqHXm4LSUBf7+1B2KjI71/IRGRTgzqWtiScxoPz8rA0dNFeGlEd4zun+jXBv9ERJ4wqP30acYhPPfZZjSuG415D12OPomNjS6JiGyKQe2j4jIn/rZ0G2atO4DLOzTFW/ekoFn9OkaXRUQ2xqD2weH8Qjw8awM2HszHQ1d1wMTrkhEVadkZjkRkEQxqndbsOo4JczJRUubCe/f2wQ09WhldEhGFCQa1FyKC6d/twbQV2egYXx/T7+uLjvH1jS6LiMIIg9qDs0Wl+OOCjVix9RiG9WyFV0f2RL06/JYRUWgxddzYeewsxs3MwP6TBfjzsC4YM7g9p94RkSEY1BqWbDyMZz/dhHp1ojD7wf7o36Gp0SURURhjUFdR6nTh5WXZ+GjNXqS2a4x37u2DFg1jjC6LiMIcg7pc7pkiPDJ7A37edwoPDEzCn27sgugoTr0jIuMxqAH8vO8kxs/agHNFZXjjrt64pbfxGz0REVUI66AWEXy8Zh+mLNuOtk3qIm1MfyS3bGB0WUREFwnboD5fXIZJizZj6cbDuLZrC7x+Ry80jHEYXRYRUQ1hGdR78s5hXFoGduWew8ShyXj4qo6IiODUOyIyp7AL6uVbjuKPCzYiOioC//l9fwy+pJnRJREReaRrWoNS6nql1A6l1C6l1KRgFxUMZU4XXlmejXFpGegYXw9LJwxmSBORJXgdUSulIgG8A+BaAIcA/KyUWiIi24JdXKCcOFeMx+ZmYs2uE7infyKev6kr6kTxFBYisgY9rY9+AHaJyB4AUErNBXALAEsE9S/HzuJ3H63HifMlePX2nrgjta3RJRER+URPUCcAOFjl/w8B6F/9QUqpsQDGAkBiYmJAiguElo1icEmLBvj30GR0T2hkdDlERD4L2NI7EZkhIqkikhofHx+op621BjEO/Of3/RjSRGRZeoI6B0DVfkGb8o8REVEI6AnqnwFcopRqr5SKBnAXgCXBLYuIiCp47VGLSJlS6lEAKwBEAvhIRLYGvTIiIgKgc8GLiCwDsCzItRARkQbu40lEZHIMaiIik2NQExGZHIOaiMjklIgE/kmVygOwP+BPXDvNABw3uogAsMt1APa5FrtcB8BrMVI7EdFcLRiUoDYjpVS6iKQaXUdt2eU6APtci12uA+C1mBVbH0REJsegJiIyuXAK6hlGFxAgdrkOwD7XYpfrAHgtphQ2PWoiIqsKpxE1EZElMaiJiEzO9kFth4N5AUAp9ZFSKlcptcXoWmpDKdVWKbVKKbVNKbVVKfW40TX5SykVo5Rar5TaWH4tLxpdU20opSKVUplKqS+MrqU2lFL7lFKblVJZSql0o+sJBFv3qMsP5t2JKgfzArjbSgfzVlBKXQngHID/iEh3o+vxl1KqFYBWIrJBKdUAQAaAERb9mSgA9UTknFLKAWA1gMdF5CeDS/OLUuopAKkAGorIcKPr8ZdSah+AVBGx0mIXj+w+oq48mFdESgBUHMxrOSLyPYCTRtdRWyJyREQ2lP/3WQDbceFcTsuRC86V/6+j/B9LjnyUUm0ADAPwgdG1UE12D2qtg3ktGQp2pJRKApACYJ3BpfitvF2QBSAXwNciYtVr+ReAZwC4DK4jEATAV0qpjPJDty3P7kFNJqWUqg9gIYAnROSM0fX4S0ScItIbF84S7aeUslxbSik1HECuiGQYXUuADBaRPgBuAPBIedvQ0uwe1DyY14TK+7kLAcwSkUVG1xMIIpIPYBWA6w0uxR+DANxc3tudC+AapVSasSX5T0Ryyv+dC+AzXGiBWprdg5oH85pM+Q24DwFsF5F/GF1PbSil4pVSceX/HYsLN62zDS3KDyIyWUTaiEgSLvwdWSkiow0uyy9KqXrlN6mhlKoH4DoAlp4pBdg8qEWkDEDFwbzbAcy36sG8Sqk5ANYCSFZKHVJKjTG6Jj8NAnAfLozassr/udHoovzUCsAqpdQmXBgUfC0ilp7aZgMtAKxWSm0EsB7AlyKy3OCaas3W0/OIiOzA1iNqIiI7YFATEZkcg5qIyOQY1EREJsegJiIyOQY1EZHJMaiJiEzu/wM8/sUVrYyJQwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "prostate = pd.read_table(\"prostate.data\")\n",
    "prostate.drop(prostate.columns[0], axis=1, inplace=True)\n",
    "\n",
    "X = prostate.drop([\"lpsa\", \"train\"], axis=1)\n",
    "y = prostate[\"lpsa\"]\n",
    "\n",
    "regressor = LinearRegression()\n",
    "\n",
    "regressor.fit(X, y)\n",
    "y_pred = regressor.predict(X)\n",
    "\n",
    "print(regressor.__dict__)\n",
    "print(y - y_pred)\n",
    "\n",
    "plt.scatter(y, y_pred)\n",
    "plt.plot([0, 5], [0, 5])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "class LogisticRegression:\n",
    "\n",
    "    def __init__(self, learning_rate=0.001, n_iters=1000):\n",
    "        self.lr = learning_rate\n",
    "        self.n_iters = n_iters\n",
    "        self.weights = None\n",
    "        self.bias = None\n",
    "\n",
    "    def fit(self, X, y):\n",
    "        n_samples, n_features = X.shape\n",
    "\n",
    "        # init parameters\n",
    "        self.weights = np.zeros(n_features)\n",
    "        self.bias = 0\n",
    "\n",
    "        # gradient descent\n",
    "        for _ in range(self.n_iters):\n",
    "            # approximate y with linear combination of weights and x, plus bias\n",
    "            linear_model = np.dot(X, self.weights) + self.bias\n",
    "            # apply sigmoid function\n",
    "            y_predicted = self._sigmoid(linear_model)\n",
    "\n",
    "            # compute gradients\n",
    "            dw = (1 / n_samples) * np.dot(X.T, (y_predicted - y))\n",
    "            db = (1 / n_samples) * np.sum(y_predicted - y)\n",
    "            # update parameters\n",
    "            self.weights -= self.lr * dw\n",
    "            self.bias -= self.lr * db\n",
    "\n",
    "    def predict(self, X):\n",
    "        linear_model = np.dot(X, self.weights) + self.bias\n",
    "        y_predicted = self._sigmoid(linear_model)\n",
    "        y_predicted_cls = [1 if i > 0.5 else 0 for i in y_predicted]\n",
    "        return np.array(y_predicted_cls)\n",
    "\n",
    "    def _sigmoid(self, x):\n",
    "        return 1 / (1 + np.exp(-x))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "LR classification perf:\n",
      " [[88  9]\n",
      " [40 16]]\n",
      "LR classification error rate:\n",
      " 0.3202614379084967\n"
     ]
    }
   ],
   "source": [
    "heart = pd.read_csv(\"SAheart.data\")\n",
    "heart.famhist.replace(to_replace=['Present', 'Absent'], value=[1, 0], inplace=True)\n",
    "heart.drop(['row.names'], axis=1, inplace=True)\n",
    "X = heart.iloc[:, :-1]\n",
    "y = heart.iloc[:, -1]\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33, random_state=42)\n",
    "\n",
    "regressor = LogisticRegression(learning_rate=0.0001, n_iters=1000)\n",
    "\n",
    "regressor.fit(X_train, y_train)\n",
    "y_pred = regressor.predict(X_test)\n",
    "perf = sklearn.metrics.confusion_matrix(y_test, y_pred)\n",
    "print(\"LR classification perf:\\n\", perf)\n",
    "\n",
    "error_rate = np.mean(y_test != y_pred)\n",
    "print(\"LR classification error rate:\\n\", error_rate)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
